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Abstract 

Acid ionization in aprotic media is studied using Molecular Dynamics tech- 
niques. In particular, models for HC1 ionization in acetonitrile and dimethyl- 
sulfoxide are investigated. The proton is treated quantum mechanically using 
Feynman path integral methods and the remaining molecules are treated clas- 
sically. Quantum effects are shown to be essential for the proper treatment 
of the ionization. The potential of mean force is computed as a function 
of the ion pair separation and the local solvent structure is examined. The 
computed dissociation constants in both solvents differ by several orders of 
magnitude which are in reasonable agreement with experimental results. Sol- 
vent separated ion pairs are found to exist in dimethylsulfoxide but not in 
acetonitrile. Dissociation mechanisms in small clusters are also investigated. 
Solvent separated ion pairs persist even in aggregates composed of rather few 
molecules, for instance, as few as thirty molecules. For smaller clusters or 
for large ion pair separations cluster finite-size effects come into play in a 
significant fashion. 
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I. INTRODUCTION 



The ionization of an acid, HA, in condensed phases 

HA^H+(S)+A-(S) (1) 

is clearly influenced by the nature of the solvent S in which the ionization takes place. A 
molecular understanding of the processes responsible for these solvent effects requires an 
analysis of the solvation forces that bring about the ionization and examination of how 
these forces depend on the properties of the solvent molecules. The analysis is complicated 
by the fact that the proton is a quantum object and thus the theoretical description must 
account for its quantum character in determining the interactions between the proton and 
the solvent. 

In this article, we study such acid ionization process in aprotic media. We have selected 
two simple prototypic molecular solvents for this study: acetonitrile (ACN) and dimethyl- 
sulfoxide (DMSO). While these two molecules have similar dipolar properties they have 
markedly different effects on the ionization process; for example, the dissociation constant 
for HC1 in DMSO is several orders of magnitude greater than that in ACN. [|TJ Furthermore, 
aprotic solvents have been chosen in order to avoid complications due to cooperative proton 
motion in the solvent (Grotthuss mechanism) [0] that may play a role when water is the 
solvent. || 

The main focus of this work is on the determination of the potential of mean force as a 
function of the separation between the H + and Cl~ ions. This allows us to characterize the 
equilibrium solvation structure; in particular, if solvent separated ion pairs (SSIP) J|] can 
exist in either solvent so that the ionization process should be represented by 

hci — > H+| |cr — ► h+(s) + cr (S) . (2) 

Here H + || Cl~ represents a metastable configuration corresponding to a local minimum in 
the mean potential where solvent molecules reside between the ions and the ions may not 
yet be considered to be completely free ions in solution. || Our calculations permit us to 
determine properties of the solvated quantum proton and to study the differences in the 
solvation properties in both solvents. 

Since the ionization process is an example of a reaction that is strongly solvent influenced, 
it is also of interest to see how the chemistry of ionization processes is changed in finite-size 
systems. Consequently, we also examine molecular clusters with linear dimensions in the 
nanometer range to study how the ionization process chang function of the cluster 

size. It is well known that chemical reactions in clusters present some unique features 
due to the lack of translational symmetry; for example, certain charge transfer reactions 
leading to neutral ion pairs take place in confined systems only above a certain threshold 
cluster size. Acid dissociative processes occurring in aqueous clusters are also known to 
exhibit similar behavior where cooperative effects gradually increase, stabilizing the reaction 
products as the number of the components in the cluster increases. || Consequently, it is 
interesting to extend bulk studies to cluster domains to unveil peculiarities of the mechanisms 
that govern dissociative processes in aggregates containing small numbers of particles. 

The organization of the paper is as follows: Section II gives details of the interaction 
potentials used in this work, along with the molecular dynamics methods used to carry out 
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the simulations of the mean potential. Bulk phase and cluster results are contained in Sees. 
Ill and IV, respectively. The conclusions of this study are given in Sec. V. 



II. SYSTEM DEFINITION 
A. Solvent potential parameters 

The ACN [CH 3 CN] and DMSO [(CH 3 ) 2 SO] solvents considered in this study were mod- 
eled as rigid molecules with the following characteristics: ACN was comprised of three sites 
representing a united atom model for the CH3 group, C and N units in the linear molecule 
with separations d between sites in a molecule fixed by constraints || at the following values: 
d(C-CH 3 )=1.46 A and d(N-CH 3 )=2.63 A. The molecular dipole moment is ^acn = 4.14D. 
The molecules interact by site-site Lennard- Jones (LJ) and Coulomb interactions. The LJ 
parameters and partial charges were taken from Edwards et al. JL0J and interactions be- 
tween identical sites are given by = 3.3, oq = 3.4 and <tch 3 — 3.6; = 50, ec = 50 
and €ch 3 = 191; qn = —0.398, qc = 0.129 and qcu 3 = 0.269. Here and below the a values 
are expressed in A, the e parameters in K and the partial charges in units of the electronic 
charge e. The usual arithmetic and geometrical means were used to determine respectively 
the a and e parameters for the cross interactions. 

DMSO was modeled as a rigid four-site molecule with the following internuclear sep- 
arations between constituents fTJ: d(S-CH 3 )=1.799 A, d(S-0)=1.48 A, d(0-CH 3 )=2.657 



A and g?(CH 3 -CH 3 )=2.685 A. The molecule is nonlinear with the following angles among 
the constituents: #(CH 3 -S-CH 3 )= 96°34', #(CH 3 -S-0)= 106°39' and the angle that the 
OS bond makes with the CH 3 -S-CH 3 plane is 64°30'. The molecular dipole moment is 
Udmso — 4.42ZJ. The LJ parameters and partial charges were taken from Rao and Singh 



I2|j and for interactions between the identical sites are given by ao = 3.3, as = 4.0 and 
a CH3 = 4.03; e = 33.2, e s = 101.7 and e C H 3 = 80.5; q = -0.459, q s = 0.139 and 
qcH 3 = 0.160. Parameters for cross interactions were determined as described above. 



B. HC1 and solvent-solute interactions 

We consider a simple model for the ionization of HC1 in the above solvents. In the gas 
phase the lowest energy potential energy surface corresponds to dissociation into neutral 
atoms while in solution, the molecule dissociates into ions due to solvation forces. A full 
quantum calculation of the potential energy surface appropriate for the ionization in solution 
must account for solvent degrees of freedom. We considered a much simpler model. The 
interaction potential energy between the H and CI atoms leading to dissociation into ions 
was determined from a density functional calculation for a Cl~ ion in the field of a point 
positive charge a fixed distance away. By appropriate distribution of the basis functions 
describing the electronic density of the complex in the neighborhood of the chlorine atom, 
the "bare" molecule was artificially constrained to dissociate into Cl~ and H + ions at large 
separations, while the charge density was allowed to redistribute at shorter separations. The 



computation was performed at the Generalized Gradient Approximation level [13| using a 
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flexible basis set for the CI . JT4J Results of these calculations as a function of the interionic 
distances r were fitted using the following potential energy function 



e 2 



V H+ ci- (r) = A exp(-Ar) - - , (3) 



r 



where e is the electron charge, A = 1.273 x 10 4 kcal/mol and A = 4.159 A -1 . Note that in this 
functional form the Coulomb interaction does not depend on the partial charges on the ions 
determined from the computed electronic density. Nevertheless, in the interactions between 
the ions and the solvent sites, the variation of the ionic charge with the distance r between the 
CI - and H + ions was taken into account. In this way, the changes in the electronic density 
of the HC1 molecule during the dissociation process was taken into account. ^From our 
density functional calculations, the charge on the Cl~ ion for r > 1 A could be represented 
approximately by the formula 

qc L {r) = -1 + 0.905 e^ r , (4) 

with 7 = 0.48 A" 1 , g#(r) = —qciij)- In performing the calculations, we have made no 
attempt to reproduce with ([3]) the actual value of dissociation energy of HC1 in vacuo or 
the exact value of the equilibrium distance for the HC1 chemical bond. In fact, our simple 
approach based on density functional theory was used only as a reasonable interpolative 
scheme between two well defined electronic structures for reactants and products states: the 
covalent character of the intramolecular bond at short interatomic distances and the ionic 
character of the dissociated species. Although this is a crude approximation, it captures 
the importante features of the solvation effects in different environments, especially at the 
larger separations of interest in this study where the species are well approximated by H + 
and Cl~ ions. 

Using the charge density extracted from this calculation, the interactions between the 
Cl~ and H + ions and the solvent molecules were taken to be site-site LJ plus Coulomb 
interactions. The LJ parameters for these ions were taken from data on Cl~ in H2O given 
in Rossky et al. [IS| using the normal geometric and arithmetic means and the known 



parameters for SPC water. The values used in the present calculation are Gci- = 3.93 
and o~h+ = 0.35; eci- = 0.830, and e#+ = 0.155 and, as earlier, the LJ parameters for 
cross interactions between these ions and the molecular solvent sites were determined from 
geometric and arithmetic means. 

C. Simulation details 

The proton was treated quantum mechanically using Feynman's path integral formula- 
tion of quantum mechanics. [TO In this representation the proton is represented by a ring 



polymer with harmonic bonds between the polymer monomers, which also interact with the 
other classical particles in the system. [0,0 The effective potential of the system is 

V eff = V d (R N ) + ^ f> - r m ) 2 + ±V p (r P , R N ) , (5) 
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where is Boltzmann's constant times temperature, m p represents the proton mass, H N 
refers to the set of classical coordinates and r p to the proton monomer coordinates with P 
the number of monomers (r p+1 = r 1 ). Here V^(R ) is the potential energy of all classical 
particles and V p (r p , H N ) is the interaction energy between the proton monomers and the 
classical particles. We have taken P = 20 in order to accurately represent the proton but 
the results were checked against simulations with P = 40 and no significant changes in the 
equilibrium averages were seen. 

Equilibrium averages were computed from time averages over a fictitious molecular dy- 
namics (MD) generated by the following Hamiltonian: 

p 1 N 1 

# = E ~ 2 rn^ + £ -MjR| + V eff , (6) 
i=i j=i 

where m* is the fictitious mass assigned to the proton monomers. We have taken m* = 30 
atomic units. Given this choice the oscillation period for the harmonic forces in the polymer 
is 

2vr / m* \ 1/2 m 

r osc = — = 2tt -— (3h. 7 
uj \ Pm p J 

Using the values given above, at T=200 K and T=293 K the period is t osc = 207 fs and 
t osc = 109 fs, respectively. 

The bulk phase molecular dynamics simulations were carried out on a system composed 
of 142 solvent molecules plus the two ions confined to a periodic box with sides L. The LJ 
interactions were cut off at half the size of the simulation box. A quartic spline interpola- 
tion was used to make the the Coulomb interactions go smoothly to zero at L/2 in a l A 
region. The box sizes were L=23.2 A for ACN and L=25.78 A for DMSO giving densities 
of 52.2 cm 3 /mol and 71.7 cm 3 /mol, respectively. The bulk simulations were carried out at 



constant temperature using Nose dynamics. [18 1 Two Nose themostats were used, one for 



the proton polymer and the other for the remaining classical particles. The temperatures of 
both thermostats were fixed at 293 K. The MD integration was performed using the Verlet 
algorithm |T^] with a time step of 2 fs. Note that with this time step there are between 50 
and 100 time steps per oscillation period of the proton monomers. 

The cluster results were obtained from time averages over constant energy molecular 
dynamics simulations. The clusters were equilibrated for 50 ps using constant temperature 
(Nose) molecular dynamics after this period the thermostat was switched off. Time averages 
were then determined from 4-5 ns constant energy MD trajectories. Clusters with sizes 
ranging from n = 4 to n = 30 solvent molecules were studied. The average temperature was 
200 ± 20 K for the small clusters with much smaller fluctuations for the larger clusters. 



III. ACID IONIZATION IN BULK SOLVENTS 

Our study of acid ionization will be primarily concerned with the change in character of 
the ionization process, as reflected in the potential of mean force or free energy for varying 
separation of the ions, as a function of the different solvent environment. In particular, we 
will be interested in describing (i) the magnitude of the resulting free energy barriers for 
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the dissociation processes in both solvents and (ii) the main features that characterize the 
solvation structure of the complex as the reactive process takes place. 



A. Potential of mean force 



The mean potential as a function of the distance between the ions was obtained from 



the mean force on the ion pair. |20| The proton position was defined as the centroid of the 
proton polymer, [pl[p2 



P 



i=l 



and thus the distance between the CI and H + ions is 



IR, 



ci 



(8) 



(9) 



with Re; the position of the CI ion. The mean force at a fixed ion pair displacement r was 



calculated as 20 



F(r) = (J2F l -r) r + 



i=i 



(3r 



(10) 



where F, is the force exerted by the Cl~ and the solvent particles upon the i-th proton 
monomer, r represents the versor along the interionic distance and (....) r represent time 
average computed using the constrained reaction coordinate dynamics ensemble. [23[] The 
mean potential W(r) was then found by integration of the mean force using 



F(r) 



dW(r) 
dr 



(11) 



Figure 1 shows the potential of mean force for the ion pair in bulk ACN and DMSO 
solvents. For reference the "bare" HC1 potential in the absence of solvent is also shown. The 
potentials have been shifted to correspond at their first minima. (We term this the contact 
ion pair state (CIP). In this case there is no distinction between the bound molecular species 
and the contact ion pair state.) We note the pronounced differences between the ACN and 
DMSO solvation. For example, at r max = 4.0 A, the maximum in the DMSO mean potential, 
the curves differ by (3AW(r max ) = W A cN{r max ) - W D MSo(r max ) = 10.9. From these curves, 
we can easily extract the ratio of dissociation constants K d 



I (1 - e DM so(r)) e-^^M 



dr 



K d (DMSO) _ [l OAcv(r) e^^M dr ] . , ■..„„,,„„. ,,, . ^ ? 

K d (ACN) ~ [J e DM so(r) e"^«°M dr] [/ (1 - Q A CN(r)) e -/roWr) dr ] 2 ^ 
where Oj(r) is unity for those r values which correspond to configurations of solvent i exhibit- 



ing associated ion pairs; for other values ©i(r) is zero. Equation (12) can be approximated 
by 



K d {DMSO) ^ .[ W e -PWACN(r) r 2 dr 
K d {ACN) ~ Jq™* e -^ w ^so(r) r 2 dr 



-/3AW(r max ) 



10 £ 



(13) 



where we have dropped terms of order one and have assumed that the boundary between 
associated and dissociated ionic states lies at the first maximum beyond the CIP well. Our 
simulations not only predict the correct experimental trend but also the result is in rough 
accord with the large experimentally observed solvent effect on the acid ionization [p]], where 
[K d (DMSO)/K d (ACN)] exp = 6.9 10 6 . Finally, we note that the individual dissociation 
constants for both solvents differ by several orders of magnitude from the corresponding 
experimental values. This is expected in view of the simplifications introduced in the design 
of the potentials. In fact, differences in a few tenths of a kcal in the potentials lead to 
considerable variations in calculated dissociation constants. Since our emphasis is on the 
study of relative, quantitative changes in the ionization process as a function of solvent type 
and environment, we have made no attempt to tune our simple potential model to describe 
all quantitative aspects of the equilibrium structure. 



B. Bulk solvation structure 

The profiles of the potential of mean force presented in Fig. 1 show that there are no 
clearly identifiable solvent separated ion pairs for ACN, as signaled by a secondary minimum 
in the mean potential. There is only a shallow minimum at r s = 5.0 A giving a barrier 
between the SSIP and CIP states of AW A cN(bulk) = W A cN(r max ) - Wacn^s) — 1.0 
kcal/mol which is within the statistical uncertainty of our calculations. In contrast, the mean 
potential for DMSO shows a fairly deep secondary minimum at r s = 7.5 A with a barrier 
separating the CIP and SSIP states of AW DM so(bulk) = W DM so{r m ax) - W DM so{r s ) = 
5.4 kcal/mol. In the DMSO calculation, the curve was computed from a constrained MD 
simulation as described above. The filled circles were determined using the formula 

P(r) = P.c-^W , (14) 

where P(r) is the probability density of r and P u is the uniform distribution. The probability 
density P(r) was estimated from a histogram of the separation between the ions obtained 
from a long unconstrained MD trajectory initiated near the mean potential secondary min- 
imum at r s = 7.5 A. Note that most of the probability density is confined between 6 — 10 
A with no escape to either contact ion pairs or free ions observed in the course of the 4 ns 
MD simulation. 

An idea of the average solvation structure in DMSO can be obtained from Fig. 2, where 
we show a contour plot of the solvent charge density, n q (p,z), for an interionic separation 
of r = 7.5 A, which corresponds to SSIP separation configuration. The charge density is 
defined as 

p B n q (p, *) = <££ ^(P? ~ P)^t - z)) r , (15) 

where p B represents the solvent bulk density and the angular brackets (• • -) r denote an 
ensemble (or time) average over configurations with the ion pair separation at a given value 
of r. In (|15|), we have selected a frame of reference centered on the ion pair, which is taken to 
lie along the z direction in a cylindrical coordinate system (p, z, 4>); (pf, zf) and q a represent 
the cylindrical coordinates and charge of site a in the i-th solvent molecule, respectively. 
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The CI - and H + ions are denoted by semicircular shaded regions in the bottom of the figure 
while the + and — signs indicate regions of high charge density of the corresponding sign. 
The picture provides information about both the locations of high solvent density as well as 
the average orientation of solvent molecules in the vicinity of the ion pair. Figure 2 clearly 
shows that the solvation structure for DMSO is energetically dominated by the electrical 
coupling between the proton, a bare charge of small size, and the negatively charged oxygen 
sites of the solvent. This is indicated by the prominent " — h +" charge density maxima to 
the right and above the proton; the molecular dipoles are oriented as might be expected, 
with the negative ends of the molecules lying near to the H + side of the ion pair complex 
and the positive ends of the other solvent molecules oriented towards the Cl~ side. However, 
the positive charge density near the Cl~ ion is more diffuse indicating weaker solvation and 
less charge localization. The analysis of the solvation structure in ACN for a similar ion pair 
separation shows essentially the same qualitative features as previously described for DMSO, 
so a comparative analysis of the charge density profiles is insufficient to understand the 
differences in ionic stabilization that these solvents present. Perhaps this is to be expected 
since the solvent contribution to the potential of mean force is the result of a subtle interplay 
between packing effects governed by short-range, repulsive forces related to the molecular 
shape, and solute-solvent dipolar forces which depend on overall distribution of molecular 
electronic density. 

The analysis of the solvent structure in the neighborhood of the transition region between 
CIP and SSIP ion pair configurations presents interesting differences for the two solvents. 
We examine the region where the steep attractive branch of the mean potential reaches its 
local maximum and levels off. Consider the dipolar density, 

Pb z) = ^—S(pi - p)S(zi - z)) r , (16) 
i 2n P 

in the pz-plane. Here fa is the dipole moment of molecule % and, as above, the average is taken 
with the distance between the ions fixed at r. Figures 3(a) and (b) show n M (p, z) for r = 2.64 
for ACN and DMSO, respectively. For ACN one observes two large peaks with opposite sign 
in the vicinity of H + . This indicates that the ACN molecules have a considerable degree of 
orientational variability so that when the molecular dipole moment is projected onto the z 
axis of the ion pair the projection can take either sign. In contrast, the dipolar density for 
DMSO in Fig. 3(b) shows no such effect implying that the DMSO molecules are much more 
rigidly ordered in the first solvation shell as indicated in the schematic representation of the 
structure shown in the bottom parts of the figures. The results suggest that in ACN the 
local solvent density fluctuations are sufficiently large to promote dipolar orientations that, 
in average, lead to a less effective dielectric reactive field, lowering the energetic cost to pass 
from SSIP to CIP states. However, in DMSO a much more structured three-dimensional 
solvent network exists that hinders molecular rotations and leads to a stronger coulombic 
coupling between the proton and the negatively charged oxygen site and a larger free energy 
barrier. 
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IV. ACID DISSOCIATION IN CLUSTERS 



Compared to bulk environments, dissociation processes in small clusters present some 
distinctive features that may lead to significant variations in the overall dissociation energy 
of the ion pair. Perhaps the simplest question to be considered is how large should a cluster 
be to exhibit bulk-like behavior? The answer will depend on the particular stage of the 
dissociation process one is interested in describing. In this study we focus our attention on 
the small-intermediate interionic distance regime; i.e, those characteristic of CIP and SSIP 
configurations seen in bulk phases. For interionic distances which are comparable to typical 
linear dimensions of the cluster, surface forces introduce new effects in the dissociation 
mechanisms. |24| 

In Fig. 4 we present the potential of mean force for HC1 dissociation in clusters of ACN 
and DMSO. We first observe that the results for the mean potential for an n = 14 ACN 
cluster, for interionic distances that are not too large, lie close to those for bulk solvents, 
indicating that even for this relatively small aggregate the solvent reactive field is comparable 
to that of bulk environments and that surface effects play only a minor role in determining 
the functional form of W(r). Very small clusters exhibit quite different behavior as is 
evident from the results for W(r) for an n = 4 ACN cluster shown in the figure. The typical 
bulk phase plateau regime at moderate interionic distances has disappeared since for such 
small aggregates there is insufficient space for independent solvation of the individual ions. 
Figure 4 also shows results for the mean potential for a DMSO cluster with n = 30 at 
temperature T = 200 K. Our results confirm that solvent separated ion pairs still persist 
in clusters of this size although the free energy barrier between the CIP and SSIP states is 
smaller, AW£>Mso( n — 30) = 1.8 kcal/mol. 

Several interesting features of the dissociation mechanisms in ACN which are common to 
bulk and clusters environments can be seen from typical snapshots of cluster configurations. 
In Fig. 5 we present two often-encountered cluster configurations for ACN clusters with 
n = 30 where the ion pair separation is fixed at the transition value corresponding to the 
sharply rising part of the mean potential discussed above. The ACN molecules lying closest 
to the H + ion are rendered in dark colors for contrast while the remaining solvent molecules 
are lightly shaded. One sees that while on average four ACN molecules strongly solvate the 
H + ion, their orientations may be such that the projection of their dipole moments on the z 
axis can be positive or negative. In contrast, DMSO clusters do not exhibit this dual dipolar 
orientation. 

Conformational equilibria between clusters of distinctive geometry occur as we approach 
the limiting regime where clusters split into two or more independent aggregates. pl| The 



interionic distance at which these phenomena are observed dependens on the cluster size. 
For the very small n = 4 ACN cluster case shown in Fig. 6 this occurs for an ion pair separa- 
tion r = 5.29 A. At this ion pair separation one sees the two types of solvent configurations; 
namely, (a) three molecules strongly solvating the H + ion with the fourth solvent molecule 
lying farther away and solvating the Cl~ ion. The other prominent configuration (b) occurs 
when all four solvent molecules strongly solvate the H + ion. This is a transitional distance: 
for smaller separations one always sees the solvation structure given by (a) while for larger 
separations one always sees that corresponding to (b). If one considers larger clusters anal- 
ogous phenomena exist, however here we have to consider larger ion pair separations. As an 
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example, Fig. 7 shows results for a n = 14 ACN cluster with the ion pair separation fixed 
at the very large value r = 27 A. The equilibrium structure is determined by configurations 
(a) where varying numbers of solvent molecules solvate the H + and Cl~ ions but these two 
solvated ions form subclusters which only weakly interact. However, there also exist equi- 
librium configurations (b) where long "strings" of solvent molecules form a link or "rubber 
band" connecting the two ions leading to a significant attractive force. This is a distinct 
cluster effect which has no counterpart in the bulk. Such structures were observed earlier for 
ions of like charge in water. [[24[] However, we remark that in the regime of ion pair separa- 
tions investigated there is no hysteresis in the mean potential indicating that the equilibrium 
configuration space is being sampled fully even for the "rubber band" configurations shown 
in Fig. 7. 

Two examples of cluster structure with the ion pair at r = 7.5 A, the SSIP configuration 
for DMSO, are shown in Fig. 8, one for ACN (top) and the other for DMSO (bottom). 
The noteworthy features of both panels are the obvious strong orientational order near H + 
for both solvent species. However note that although there is little discernible difference 
between the two solvent cases, SSIPs exist for clusters of DMSO but do not for ACN. These 
two pictures further illustrate the fact that the cluster shape fluctuations are strong and 
these fluctuations must be taken into account in any model of the cluster solvation dynamics. 
Furthermore, one sees that surface forces do manifest themselves in the orientational order of 
the solvent molecules on the cluster surface which tend to lie parallel to the surface. However, 
as observed earlier, even for fairly small clusters (n = 14) there are only rather small 
quantitative effects on the mean potential and such surface force effects become important 
only at longer ion pair separations or for smaller cluster sizes. 

V. CONCLUSIONS 

The calculations of the ionization of HC1 in ACN and DMSO solvents indicated some 
of the solvation features that make the ionization process rather different in these solvents 
in spite of their similar dipole moments. The differences could be ascribed to the details 
of the charge distributions in the molecules and their effects on the quantitative aspects 
of the solvation structure. The gross features of the solvation structure were shown to be 
similar for both solvents. However, in spite of the apparent similarities in the local solvation 
structure, DMSO supports solvent separated ion pairs while the tendency to form such pairs 
is considerably reduced or absent in ACN. 

The cluster environment has an important influence on the ionization process, especially 
for clusters with fewer than n = 14 molecules. For ion pair separations as large as that for 
the SSIP state, the cluster potential of mean force closely corresponds to that of the bulk, 
even for clusters with as few as n = 14 solvent molecules. For much smaller clusters, as 
might be expected, there are significant changes in the mean potential reflecting the finite 
size of the cluster environment. Of course, as the ion pair is stretched the finite size of the 
cluster gives rise to much more dramatic effects on the mean potential as reflected in the 
"rubber band" -like solvent configurations that lead to attractive forces between the ions that 
cannot exist in the bulk. Since these forces correspond to an equilibrium cluster and involve 
averages over all cluster configurations with the ion pair separation fixed, they may have 
little relevance for ion pair dissociation processes in clusters since the time scale to establish 
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such equilibrium structures may be long compared to typical times for dissociation. It would 
be interesting to investigate the dynamics of dissociation in clusters in this context. 

In order to estimate the importance of the explicit incorporation of the quantum nature 
of the proton in our calculations, we performed a series of simulations in an n = 14 ACN 
cluster where the proton was treated as a classical point charge. The results for the computed 
potential of mean force are included in Fig. 4; they clearly show the necessity to treat the 
proton quantum mechanically in order to accurately represent the activation free energy 
for ionization in these solvents. Classical mechanics underestimates the barrier for passage 
from the CIP to SSIP or free ion states. In this respect, the present behavior is opposite to 
that found in proton transfer reactions where quantum tunneling reduces the height of the 
effective free energy barrier. [2(| Quantum dispersion tends to delocalize the proton charge 
leading to a less effective ionic solvation, as well as smaller solvation energies for dissociated 
ions in comparison to results from classical models where the proton is taken to be a localized 
point charge. 

The simulations presented in this paper show that acid ionization in molecular clusters 
can mimic that in the bulk even for quite small clusters. However, for clusters below a 
certain size, about n = 14 in this study, distinctive finite-size effects influence the ionization. 
The study points to some of the new chemistry that can occur as a result of environmental 
changes on the microscopic and mesoscopic levels. 
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FIGURES 




FIG. 1. Potential of mean force in units of kcal/mol for bulk phase systems: solid line, bare 
potential; dashed lines, ACN; dashed-dotted lines, DMSO; open circles correspond to unconstrained 
trajectories results for DMSO and ACN at small distances; filled circles correspond to unconstrained 
trajectories results for DMSO at large distances. The distance r is in units of A. 
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FIG. 2. Contour plots of the charge density n q in the pz-plane with the ion pair in the SSIP 
configuration at r = 7.5 A. Top and bottom panels correspond to DMSO and ACN, respectively. 
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FIG. 3. Dipolar density (written as fi on the ordinate) with the ion pair at r = 2.64 A: (a) 
ACN and (b) DMSO. Schematic representations of the local solvent structure are shown in the 
botton parts of the figures. 
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FIG. 4. Potential of mean force for clusters: solid line, bare potential; dashed line, ACN cluster 
n = 14; dotted line, ACN cluster n = 4; dotted-dashed line, DMSO cluster n = 30; short-dashed 
lines, ACN cluster n = 14 treating the proton classically. Units are the same as in Fig. 1. 
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FIG. 5. Two configurations of an n = 30 ACN cluster with r = 2.64 A. The molecules strongly 
solvating H + are heavily shaded and the two configuration illustrate the different orientations of 
the solvent molecules around H + . 
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